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Introduction 


The  location  of  buried  tunnels  is  part  of  an  ongoing  pro- 
gram at  DARPA  with  many  different  disciplines  and  techniques 
involved  in  the  effort.  Geo-Centers,  Inc.  effort  to  locate 
buried  tunnels  was  originally  funded  by  DARPA  for  a period  of 
I five  months  beginning  in  June  1976.  The  major  thrust  of  that 

j effort  was  the  theoretical  investigation  of  special  computer 

assisted  detection  methods  to  locate  tunnels.  The  basic  elec- 
; tromagnetic  signal  was  theoretically  generated  using  the  prin- 

ciples of  video  pulse  radar.  This  radar  technique  uses  a short 
time  limited  electromagnetic  pulse  radiated  by  a movable  an- 
tenna platform  at  the  earth-air  interface.  The  objectives  of 
the  original  tunnel  detection  study  were  successfully  realized; 
namely,  the  development  of  a rapid  and  sensitive  signal  pro- 
cessing technique  which  was  tested  using  a number  of  realistic 
models . 

► 

I 

The  major  emphasis  of  the  Geo-Centers'  program  was  the 
development  of  a signal  processing  method  to  detect  buried 
tunnels  using  video  pulse  radar  probes.  Video  pulse  radar  is 
used  to  transmit  and  receive  a short  unidirectional  electro- 
magnetic pulse.  The  return  pulse  in  turn  produces  a two  dimen- 
‘ sional  sounding  or  picture.  Theoretical  soundings  were  generated 

based  on  a simple  scattering  model  relevant  to  long  cylindrical 

tunnels  that  are  buried  in  granite.  A fast  and  sensitive  pic- 
I 

I ture  processing  algorithm  was  evolved  for  the  specific  problem 
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of  subsurface  tunnel  detection . The  utility  of  the  algorithm 
has  been  tested  in  a variety  of  situations  including;  (1)  mul- 
tiple tunnels,  (2)  large  noise  or  clutter  background  values, 
and  (3)  interfering  or  blocking  layers.  In  all  cases  the  detec- 
tion algorithm  is  unusually  stable  and  dependable.  The  present 
method  requires  a human  operator  for  final  recognition. 

The  basic  theory  for  generating  the  subsurface  soundings 
is  related  to  the  "geometrical  optics"  approach  of  diffraction 
theory.  The  transmission  and  reflection  coefficients  at  the 
tunnel  walls  are  given  approximately  by  the  plane  wave  values. 
Geometrical  ray  paths  are  used  to  compute  the  tim.e  for  each 
tunnel  return.  This  method  of  generating  subsurface  soundings 
is  quite  adequate  except  in  the  lower  frequency  range  of  the 
incident  pulse.  However,  it  is  the  high  frequency  part  of  the 
pulse  that  is  concerned  with  resolving  tunnel  location. 

The  basic  tunnel  detection  scheme  involves  the  application 
of  a new  concept  in  signal  detection.  In  this  scheme,  a "geo- 
metrical filter"  has  been  used  to  detect  the  special  tunnel 
returns  that  are  found  on  the  two  dimensional  soundings.  The 
actual  method  used  to  construct  the  filter  is  based  on  real 
space  reconstruction  methods  used  to  generate  three  dimensional 
pictures  from  a series  of  two-dimensional  projections.  One- 
dimensional projections  are  used  to  reconstruct  the  sounding. 
They  are  chosen  so  as  to  give  superior  tunnel  representations. 
The  reconstruction  algorithms  act  as  a shape  sensitive  filter 
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to  emphasize  only  the  tunnel  attributes  of  the  soundings  and 
are  insensitive  to  other  objects  and  clutter  on  the  original 
sounding.  Once  the  soundings  have  been  processed  by  the  "geo- 
metrical filter,"  the  tunnel  locations  are  clearly  evident  to 
the  observer. 

A.  The  Application  of  Pattern  Recognition  Techniques  for 
Tunnel  Location* 

In  our  recent  study,  representative  subsurface  soundings 
were  generated  using  a modified  geometrical  optics  method  for 
scattering  by  dielectric  bodies.^  The  basic  building  block 
used  in  this  theory  is  the  time  domain  pulse  reflected  from  a 
semi-infinite  medium.  Pulses  scattered  from  a cylinder  are 
decomposed  into  a superposition  of  reflections  from  the  sur- 
faces of  a buried  cylinder  used  to  model  the  tunnel.  Geometri- 
cal paths  were  used  to  compute  the  times  for  the  individual 
reflections.  Steady  state  transfer  functions  were  superimposed 
through  the  Fourier  integral  to  produce  the  time  domain  signals. 
A sample  steady  state  transfer  function  is  shown  in  Figure  1 
and  a pulse  reflected  from  a tunnel  buried  in  granite  is  given 
in  Figure  2.  A Gaussian  incident  pulse  of  specified  half-width 
is  assumed  to  be  propagated  into  the  earth  at  ground  zero. 

Values  of  the  electrical  parameters  were:  conductivity  o = 10  ^ 

and  10  and  relative  dielectric  constant  ~ ^ ^ two 

dimensional  sounding  formed  by  moving  a probing  antenna  over 


Material  presented  in  this  section  is  considered  pro- 
prietary by  DARPA  and  Geo-Centers,  Inc. 
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the  earth  near  a buried  tunnel  is  shouT.  in  Figure  3.  This 
figure  can  be  understood  by  referring  to  the  actual  physical 
arrangement  shown  in  Figure  4.  V.’hen  the  antenna  is  on  the  left 
the  physical  distan.ce  to  the  tunnel  is  greater  than  when  the 
sound  is  overhead.  The  high  numbers  in  Figure  3 represent 
tunnel  reflections  and  present  a generally  hyperbolic  trace. 

The  small  sketch  on  Figure  4 shows  this  diagrammatically . Using 
the  maximum  nuir.ber  that  appears  on  Figure  3 a gray  scale  from 
zero  to  ma.ximum  optical  density  can  be  constructed.  This  is 
also  illustrated  on  the  right  in  Figure  3.  Note  that  the  den- 
sity variations  correspond  to  the  maximum,  and  minimiUm.  amplitudes 
on  the  time  response  curve  shown  in  Figure  2.  More  detail  than 
shown  on  Figure  3 is  possible  if  a gray  scale  with  narrower  in- 
dividual limits  is  chosen.  Although  it  is  very  easy  to  recog- 
nize a tunnel  from  this  picture,  noise,  clutter,  blocking  layers, 
and  multiple  targets  m.ake  practical  identification  very  diffi- 
cult. Conventional  methods  for  tunnel  detection  would  use  indi- 
vidual signals  and  frequency  sensitive  filters. 

The  basic  tunnel  detection  scheme  used  in  our  study  involved 
a new  concept  in  interpreting  soundings.  A "geometrical  filter" 
was  used  to  detect  tunnels  from  the  sounding.  Inverse  recon- 
struction might  also  be  a good  name  for  the  method.  Real  space 
2 

reconstructions  have  been  used  recently  to  reconstruct  a three 
dimensional  im.age  from  a set  of  two  dimensional  projections. 

As  with  classic  Fourier  reconstructions,  the  more  views  of  an 
object  that  are  available,  the  greater  the  definition  of  the 

GEO-CENTERS,  INC. 


1 

I 


21 


MODEL:  SINGLE  TUNNEL 
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reconstruction.  Tunnel  identification  is  achieved  as  follows: 

1.  Two  pseudo-projections  (each  one  dimensional)  are 
formed  at  two  different  viewing  angles  from  the 
sounding . 

2.  The  pseudo-projections  are  used  to  reconstruct  the 
sounding . 

The  angles  for  the  pseudo-projections  are  chosen  so  as  to 
reconstruct  the  tunnel  aspects  of  the  sounding  to  a high  degree 
of  accuracy.  Intervening  layers  and  clutter  are  reconstructed 
to  a low  degree  of  accuracy  with  the  reconstruction  method  (or 
algorithm) . The  result  is  that  a "geometrical"  filter  has  been 
developed  which  basically  sees  only  tunnel  tracings. 

The  final  reconstructions  show  the  details  of  minor  pulse 
pealcs  which  are  related  to  the  tunnel  size  and  the  surrounding 
material.  The  tunnel  has  a unique  pattern  which  is  emphasized 
by  the  reconstruction.  Note  that  the  sounding  picture  recon- 
struction is  very  successful  at  retaining  tunnel  attributes  in 
the  presence  of  large  amounts  of  extraneous  signals. 

Figures  5,  6 and  7 are  geometrical  filter  versions  of  the 
single  tunnel  sounding  of  Figure  3.  The  pseudo-projection 
angles  are  ±14®,  ±20°,  ±30°.  In  all  cases  the  numerical  (left) 
and  grey  scale  (rignt)  outputs  show  two  dominant  features;  a 
central  core  plus  two  oblique  skirts.  The  grey  scale  dots  are 
actually  maxima  or  minima  since  only  the  absolute  value  of  the 
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MODEL:  SINGLE  TUNNEL 

PHI  = 30' 


Cumulative  Output  Grey  Scale  Output 

Figure  7 Subsurface  Sounding  after  Geometrical  Filtering  (PHI  = 30°) 


signal  is  involved.  The  most  useful  output  is  Figure  6 since 
the  detailed  attributes  of  the  skirt  are  most  clearly  shown. 
These  skirts  contain  the  basic  information  about  the  size  and 


14 


shape  of  the  tunnel.  The  central  core  region  gives  basic  tar- 
get information  about  the  physical  location  of  the  tunnel.  The 
computer  grey  scale  sounding  for  a single  tunnel  with  a layer 
and  lower  blotch  is  shown  in  Figure  8.  A geometrically  fil- 
tered output  for  this  case  is  shown  in  Figure  9.  Note  that  only 
the  tunnel  attributes  are  preserved  by  the  picture  detection 
process . 

A noise  equal  to  25  and  50%  of  the  signal  tunnel  signal 
was  added  to  the  sounding  of  Figure  3.  The  geometrically  fil- 
tered outputs  of  Figures  10  and  11  show  that  even  in  these  ex- 
treme cases  the  pictorial  tunnel  attributes  are  evident.  Even 
when  the  added  noise  is  raised  to  100  percent  the  geometrically 
filtered  output  of  Figure  12  clearly  shows  the  power  of  the 
detection  method. 

Two  separated  tunnels  were  also  used  to  test  the  detection 
algorithm.  The  geometrically  filtered  outputs  are  shown  in 
Figures  13,  14  and  15  for  the  original  noiseless  sounding  and 
for  50  and  100  percent  additive  noise.  Again  the  tunnels  are 
clearly  evident. 

It  can  be  concluded  that  the  pictorial  detection  schemes 
is  uniquely  able  to  locate  and  classify  buried  tunnels  under 
a variety  of  conditions  including  high  noise  and  clutter  values. 
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Figure  8 Video  Pulse  Sounding  from  a Buried  Tunnel  with 
Layer  and  Blotch 
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Figure  W Geometrically  Filtered  Output  of  Two  Buried  Tunnels  Plus  Random  Noise 
(MX  OF  Fiax  Signal) 
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II . Additional  Technical  Results 

The  time  domain  signal  reflected  from  the  tunnel  was  com- 
puted from  the  steady  state  transfer  function  and  the  frequency 
spectrum  of  the  pulse.  The  transfer  functions  given  in  Figures 
1 to  8 give  the  plane  wave  reflected  signal  for  tunnels  having 

diameters  from  1 to  5 meters.  The  conductivities  range  from 
-3  -4 

1.2  X 10  to  10  S/m  and  the  relative  dielectric  constant  from 
= 7 to  = 30.  The  tunnel  is  assumed  to  be  filled  with  free 
space.  In  most  cases  the  rapidly  oscillating  high  frequency 
range  has  been  truncated. 

The  geometrical  theory  of  diffraction  used  in  the  present 
theory  involves  the  superposition  of  the  reflections  from  the 
top  granite  tunnel  interface  and  the  reflections  from  the  tunnel 
bottom.  Some  background  results  are  shown  in  Figures  9-11.  In 
Figure  9 and  10  the  pulse-  begins  in  air  and  is  traced  through 
the  granite  overburden.  Figure  11  assumes  that  the  pulse  begins 
in  granite  and  is  followed  into  the  granite  at  different  depths. 
Figures  12  and  13  show  the  initial  tunnel  reflection  from  the 
upper  tunnel  interface.  Figures  14  and  15  give  the  pulses  from 
the  lower  tunnel  interface.  Since  some  tunnel  signal  processing 
may  surpress  the  top  tunnel  signal  Figure  16  gives  the  tunnel 
return  signal  with  and  without  the  first  reflection.  Figures 
17-19  give  the  full  tunnel  returned  signals  for  2 and  7.5  meter 
tunnels  at  different  depths. 
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A.  Transfer  Functions 
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Figure  5 Steady  State  Transfer  Function,  diam  =2,  a = 1.2  x 10  S/rn, 


Ji.}  • 
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B.  Time  Domain  Pulses  (t  = pulse  width) 
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iqure  9 Gaussian  Pulse  from  Air  to  Granite,  t = 10  sec.,  o = 10  S/m 


Figure  10  Gaussian  Pulse  from  Air  to  Granite,  t = 10  sec.,  o - 10  S/m,  D - 0,i00 
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.001 


TIME  X 10 

nnel  for  Different  Diameters 


.0001 


from  Bottom  of  Tunnel,  t = 10  sec.,  o = 10  S/m 
ni . , diam  = 2 m. 


TIME  X 10 

Second  Tunnel  Pulse  from  Bottom  of  Tunnel 
e = 7,  D = 0,25,100  m. , diam.  = 2 m. 


WITH  FIRST  PULSE 


Pulse  Returned  from  Tunnel  V.’ith  and  VJithout  First  Reflection,  diam  = 2 m. 
D = 100  m.,  0 = 10““^  JjA',  £-  = 7,  t = 10~®  sec. 


.0001 


TIME  X 10 

Figure  17  Full  Tunnel  Signal,  diam  = 2 m. , o = 10“^  S/m,  e=7,  D=50,75m.,  t=10®  see. 
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Dio  = 7.5 

o"  =.001 


TIME  X 10 

Full  Tunnel  Signal,  diam  = 7.5  m.,  o = 10  ^ S/m, 
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III . Reconstruction  of  Arrays 

The  two  programs  v;hich  employed  the  reconstruction  technique 
were  BG;1  and  BJS;1.  The  first  BG;1  reconstructs  an  array  from 
two  projections,  one  vertical  and  the  other  at  an  angle  ip  from 
the  vertical.  The  other  BJS;1  reconstructs  from  two  projections 
at  ip  from  the  vertical. 

The  first  step  in  the  reconstruction  transforms  the  array 
(rectangular)  into  an  area  more  convenient  for  the  computation 
of  the  projections.  In  BG;1  the  original  rectangular  array  is 
transformed  into  a parallelogram  (Figure  1).  In  BJS;1  the  array 
is  transformed  into  an  array  whose  width  is  the  same  as  the 
original  array  but  whose  length  is  chosen  to  fit  as  closely  as 
possible  to  the  original  length.  Every  other  point  in  this 
array  is  not  part  of  the  grid  (Figure  2) . 

The  positions  of  the  new  points  are  computed  and  corres- 
ponding densities  are  found  by  using  a Cardinal  Series  Inter- 
polation technique.  That  is,  the  contribution  to  the  new  den- 
sity  from  an  old  density  point  is  given  by  — 

where  a and  b are  the  x and  y distances  between  the  points. 

The  sum  is  over  the  nearest  16  points. 

The  appropriate  projections  (P,Q)  from  these  new  densities 
are  required.  The  first  order  approximation  to  the  reconstruc- 
ted array  is  set  out  in  such  a way  that  the  Q projections  are 
satisfied.  New  P projections  are  computed  and  all  points  in  a 
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Rectangular  Array  Transformed  Figure  2 Transformation  of  Array  Points 
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given  P row  are  scaled  by  a constant  factor  so  as  to  satisfy 
the  original  P projections.  Then  the  same  is  done  for  the  Q 
projections.  This  iterative  process  is  applied  until  both  pro- 
jections are  satisfied  to  the  desired  degree  of  accuracy.  This 
usually  takes  fewer  than  four  sets  of  iterations. 


In  BJS;1  the  possible  tunnels  are  located  by  first  smooth- 
ing both  P and  Q projections  (P ( J)  = (P(J-l)  + 4P(J)  + P(J+l)/6) 
and  then  looking  at  the  intersections  of  the  maxima  of  P and  Q. 

Where  a point  determined  by  the  intersection  of  two  of  these 
projections  is  found  to  be  within  the  array,  an  average  over 
the  nearest  four  points  is  taken  and  output  along  with  the  coor- 
dinates of  the  possible  tunnel.  These  sums  are  an  indication 
of  the  "strength"  of  the  tunnel.  Figure  3 shows  the  array 
geometry  and  Figure  4 the  projection  sums. 

The  width  of  the  projection  maxima  are  an  indication  of  the 
size  of  the  tunnel,  and  finer  detail  in  the  "skirts"  of  the 

array  will  help  distinguish  tunnels  from  other  possible  targets.  > 

- 

IV.  Reflection  and  Propagation  of  Pulses  ^ 

/ 

1 

Several  programs  were  v/ritten  to  study  the  propagation  of  i 

electromagnetic  pulses  through  granite  (and  other  materials)  ! 

i 

and  the  reflection  and  transmission  of  these  pulses  from  various 
granite-air  interfaces.  The  final  program  TT.F0R;2  gives  the 
electric  field  intensity  at  the  surface  (as  a function  of  time) 
from  a Gaussian  wave  which  travels  a distance  D through  granite. 
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is  then  reflected  (three  times)  from  an  air  filled  cylinder  and 
finally  propagates  back  to  the  surface. 

The  computations  were  done  in  the  frequency  domain.  A 
steady  state  transfer  function  was  computed  and  Fourier  trans- 
formed back  to  give  a real  time  domain  signal.  The  integral 
involved  required  from  100-400  iterations  to  converge.  The 

4 

limits  of  the  integral  are  from  oj  = 10  to  co  = 2.6/Tl.  For 
lower  frequencies  the  transfer  function  was  nearly  constant  so 
the  contribution  to  the  integral  was  neglible.  For  frequen- 
cies higher  than  2.6/TT  the  Gaussian  envelope  sends  the  inte- 
gral to  zero.  The  co  integral  used  unequal  spacing  since  the 
transfer  functions  decreased  exponential  with  (o.  This  gave 

nearby  equal  intervals  areas  for  each  Ato.  (i.e.,  as  many 

7 8 6 

points  were  used  between  (o  = 10  and  10  as  between  co  =10 
and  10^,  etc.) 

Another  program  computes  the  transmitted  pulse  (air-<-granite) 
to  a depth  in  granite  to  compare  with  work  done  by  others. 

The  resultant  pulse  agreed  exactly  with  the  work  done  by  King 
and  Harrison. 

Programs  which  output  the  transfer  functions  (magnitude 
vs  frequency)  were  also  written  to  investigate  the  relation- 
ship between  tunnel  size,  pulse  width  and  particular  electrical 
characteristics  of  the  granite. 
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The  final  program  Gl.FOR;!  created  an  array  of  "density" 
points  resembling  a real  electromagnetic  subsurface  profile. 
The  return  pulses  to  a point  sweeping  across  the  earths  sur- 
face (perpendicular  to  the  tunnel)  were  transformed  to  corres- 
ponding depth  positions  and  output  on  a grid.  For  most  of 
these  investigations  the  following  parameters  were  used: 

= 7 , a = 10  ^ or  10  ^ , pulse  width  10  10  or  10 


V.  Program  Descriptions 


BJS;1  The  reflected  pulse  array  (9,43)  is  read  from 

disk,  (FOR20.DAT,  the  latest  version  is  read). 
The  computer  then  asks  for  a projection  angle 
PHI (12).  It  then  creates  an  array  whose  size 
most  nearly  fits  the  original  array  such  that 
the  projections  from  points  A(N,M)  to  A(N+1, 

M+1)  or  A(N-1,  M-1)  are  at  the  angle  ±PHI. 

New  densities  at  these  points  are  formed  from 
the  old  densities  by  using  a Cardinal  Series 
Interpolation  formula  for  the  sixteen  points 
nearest  the  new  density  point.  Projections 
(P,Q)  at  an  angle  ip  from  the  vertical  are  then 
computed  and  a zeroth  order  approximation  satis- 
fying the  Q projections.  P and  Q projections 
are  then  alternately  satisfied  until  the  array 
converges . 

The  program  outputs  the  projections,  the  recon- 
structed array,  possible  tunnel  positions  and 
then  a gray  scale  equivalent  of  the  array.  The 
tunnels  are  found  by  first  smoothing  the  P and 
Q projections  and  then  looking  at  the  intersec- 
tions of  their  respective  maxima,  1. 

After  the  gray  scale  is  output  the  machine 
waits  for  a -1  or  1 input  (12).  A -1  ends  the 
program  and  a 1 will  make  the  program  ask  for 
a new  input  angle. 

BJS;2  OUTPUTS  A GRAY  SCALE  WITH  10  shades. 
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BI;1 


BG;1 


BE;  1 


Gl.FOR; 1 


G2.FOR; 1 


G.FOR;l 


Same  as  BJS;1  but  does  not  involve  reading  the 
array  from  the  disk,  or  output  the  tunnel  loca- 
tion. It  asks  for  the  I and  J dimensions  and 
then  for  the  array.  I and  J must  be  odd. 

Reads  a (9  x 30)  array  from  the  disk  (FOR20.DAT), 
transforms  the  array  to  an  array  in  the  shape  of 
a parallelogram  and  then  proceeds  like  BJS;1 
except  the  two  projections  are  vertical  a.nd  at 
an  angle  ip.  The  output  is  the  projections,  den- 
sities and  gray  scale. 

After  the  projections  are  output  the  machine 
waits  for  an  IZ  input.  If  IZ=0  the  iteration 
is  done  two  more  times.  If  IZ=1  is  outputs  the 
gray  scale  and  the  array,  if  IZ=2  it  goes  back 
for  a new  angle  and  if  IZ=3  it  ends. 

Same  as  BG  but  does  not  read  from  the  disk  and 
also  outputs  the  original  new  densities 

An  array  is  outputed  which  corresponds  to  a scan 
across  the  earth's  surface  from  -60  to  60  meters. 
The  output  is  the  array  positioned  relative  to 
the  earth  as  it  would  appear  on  an  ESP.  The 
seventh  output  now  corresponds  to  the  position 
of  the  return  peak  from  the  x=0  reflected  pulse. 
Each  subsequent  row  is  at  a time  increment  of 
half  the  pulse  width. 

The  input  format  is : 

DEPTH  - distance  to  the  top  of  the  cylinder,  F5.1 
ITERATIONS  - number  of  iterations  in  the  inte- 
gral, 13 

SIGMA  - effective  conductivity,  F6.5 
TIME  - pulse  width,  E5.0 
DIAM  - diameter  of  the  cylinder,  II 
SCALE  - scale  factor  for  output,  12 
IW  - switch  to  include  first  pulse  IW=1  or  leave 
it  out  IW=0,  II 

Same  as  Gl.FOR;!  but  X goes  from  -40  to  40  and 
there  is  an  option  to  also  output  (IW=1)  the 
combined  signal  from  two  other  tunnels  located 
at  X = -20  and  +20  with  diam.eters  of  0 and  6. 
(This  does  not  work  correctly  because  of  how 
TT  is  used).  Note:  IW  is  not  the  same  as  in 

Gl.FOR; 1. 

Outputs  array  grid  without  the  earth  spacing  as 
in  Gl.FOR;!. 
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TT,F0R;2 


T. F0R;5 

TT.FOR; 1 
Tl.FOR;2 

TS.FOR; 1 
TS.FORj  2 

U. FOR; 1 

AC;  1 

BC;  1 
H.FOR;  1 


The  time  domain  reflected  pulse  from  a cylinder 
embedded  in  granite.  The  seventh  output  point 
corresponds  to  the  peak  of  the  pulse.  The  in- 
tegral for  each  value  of  the  time  is  performed 
by  taking  increments  of  w which  are  exponential. 

The  input  format  is: 

DEPTH  - 13  - depth  to  top  of  cylinder 
ITERATIONS  - 13  - number  of  iterations  in  the 
integral  (greater  than  100) 

SIGMA  - F6.5  - effective  condictivity 

TIME  - E5.0  - pulse  width 

DIAI4ETER  - II  - diameter  of  cylinder 

Same  as  TT.FOR; 2 but  computes  the  integral  from 
w = 10^  instead  of  lO'^  and  uses  1/2  the  value 
of  the  first  point.  (Output  varies  considerably 
from  TT.FOR; 2) . 

Same  as  TT.FOR; 2 but  without  first  pulse. 

Same  as  TT.FOR; 2 but  outputs  100  points  from 
T=0  to  T=200T1 

TT.FOR; 2 but  outputs  the  value  of  the  integral 
in  the  summing  loop  at  T = the  center  of  the 
pulse  position. 

TS.FOR;!  but  outputs  the  function  for  w = 12  x 
10®  at  T = center  of  pulse.  For  T different 
from  this,  the  output  is  markedly  different. 

Similar  to  TT.FOR;!  but  computes  the  pulse  as 
the  sum  of  three  separate  transfer  functions. 
Also,  the  cylinder  diameter  is  adjusted  to  the 
nearest  Tl/2  interval  of  delay  time.  (Phase 
difference,  in  time,  for  each  of  the  three 
pulses  is  always  a multiple  of  Tl/2) . 

Reconstruction  of  an  odd  square  array  for  an 
angle  from  the  vertical  and  the  vertical. 

Reconstruction  on  an  NxM  odd  array  as  in  ACl;. 

Computes  the  transfer  function  for  reflection 
from  a cylinder  using  first  three  reflected 
pulses.  W goes  from  10®  to  10^0,  15  points  are 
output . 

Input  format  is: 

DISTANCE  - depth  to  cylinder,  F5.1 
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HH.FOR;! 
HH.FOR; 2 


RANI ; 1 


GRl ; 1 


SIGMA  - effective  conductivity,  F6.5 
DIAMETER  - diameter  of  cylinder,  12 

8 9 

Same  as  H.F0R;1  but  W goes  from  2.6  x 10  to  10 
with  55  points  output. 

9 10 

Same  as  H.F0R;1,  but  W goes  from  10  to  10 
with  35  points  output. 

Reads  a 9 x 43  array  FOR22.DAT  adds  a random  num- 
ber to  each  element  and  writes  the  result  on 
FOR20.DAT.  The  average  random  number  is  half 
the  maximum. 

Input : 

MAX  RAN  NUMBER  - 13 

Reads  FOR20.DAT  (9  x 43)  and  outputs  gray  scale 
on  TTY.  6 shades  are  used,  0,1, 2, 3, 4, 5. 


T.FOR; 1 


S.FOR; 1 


SI. FOR; 1 


Outputs  cylinder  reflection  from  first  two  re- 
flected pulses  using  two  transfer  functions. 
The  two  pulses  are  output  separately,  not  the 
sum.  The  cylinder  has  zero  diameter. 

Reflected  pulse  going  from  granite  to  air, 
through  D meters  of  granite. 

Output  50  points  for  S.FOR; 1;  T = 4T2  200  T1 


R.FOR; 1 


Transmitted  pulse  from  air-granite  to  depth  D. 


X.FOR;  1 
Q.FOR;l 


Transfer  function  for  free  pulse  in  granite 
multi-optical  by  Gaussian  envelope  exp  (-a)2t2/2)  . 

-4 

Transmission  of  free  pulse  in  granite  a = 10  , 

t = 10“^  sec. 


P.FOR;l 

C.FOR  Plots  reflection  from  sphere. 


A.  Data  Files 


FOR21.DAT;! 

9 

X 

3 

array  < 

of 

one 

pulse  from  -40  to  40  meters. 

F0R21.DAT; 3 

9 

X 

43 

array 

of 

two 

pulses  from  -40  to  40 

meters . 

F0R21.DAT;  2 

9 

X 

43 

array 

of 

one 

pulse . 

FOR21.DAT; 4 

9 

X 

43 

array 

of 

one 

pulse  with  blotch  and 

stripe . 
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